home *** CD-ROM | disk | FTP | other *** search
/ Info-Mac 4 / Info_Mac IV CD-ROM (Pacific HiTech Inc.)(August 1994).iso / Development / Source / POVSRC / SOURCE / BLOB.C < prev    next >
Text File  |  1993-08-10  |  18KB  |  658 lines

  1. /****************************************************************************
  2. *                blob.c
  3. *
  4. *  This module contains the code for the blob shape.
  5. *
  6. *  This file was written by Alexander Enzmann.    He wrote the code for
  7. *  blobs and generously provided us these enhancements.
  8. *
  9. *  from Persistence of Vision Raytracer
  10. *  Copyright 1993 Persistence of Vision Team
  11. *---------------------------------------------------------------------------
  12. *  NOTICE: This source code file is provided so that users may experiment
  13. *  with enhancements to POV-Ray and to port the software to platforms other 
  14. *  than those supported by the POV-Ray Team.  There are strict rules under
  15. *  which you are permitted to use this file.  The rules are in the file
  16. *  named POVLEGAL.DOC which should be distributed with this file. If 
  17. *  POVLEGAL.DOC is not available or for more info please contact the POV-Ray
  18. *  Team Coordinator by leaving a message in CompuServe's Graphics Developer's
  19. *  Forum.  The latest version of POV-Ray may be found there as well.
  20. *
  21. * This program is based on the popular DKB raytracer version 2.12.
  22. * DKBTrace was originally written by David K. Buck.
  23. * DKBTrace Ver 2.0-2.12 were written by David K. Buck & Aaron A. Collins.
  24. *
  25. *****************************************************************************/
  26.  
  27. #include "frame.h"
  28. #include "vector.h"
  29. #include "povproto.h"
  30.  
  31. #ifndef min
  32. #define min(x,y) ((x)<(y)?(x):(y))
  33. #endif
  34. #ifndef max
  35. #define max(x,y) ((x)>(y)?(x):(y))
  36. #endif
  37.  
  38. METHODS Blob_Methods =
  39.   { 
  40.   All_Blob_Intersections,
  41.   Inside_Blob, Blob_Normal,
  42.   Copy_Blob,
  43.   Translate_Blob, Rotate_Blob, Scale_Blob, Transform_Blob,
  44.   Invert_Blob, Destroy_Blob
  45. };
  46.  
  47. extern long Ray_Blob_Tests, Ray_Blob_Tests_Succeeded;
  48. extern int Shadow_Test_Flag;
  49.  
  50. static int determine_influences PARAMS((VECTOR *P, VECTOR *D, BLOB *blob, DBL mindist));
  51. static DBL calculate_field_value PARAMS ((OBJECT *obj, VECTOR *Pos));
  52.  
  53. #ifndef Blob_Tolerance 
  54. #define Blob_Tolerance 1.0e-3
  55. #endif
  56.  
  57. #define COEFF_LIMIT 1.0e-20
  58. #define INSIDE_TOLERANCE 1.0e-6
  59.  
  60. /* Starting with the density function: (1-r^2)^2, we have a field
  61.    that varies in strength from 1 at r = 0 to 0 at r = 1.  By
  62.    substituting r/rad for r, we can adjust the range of influence
  63.    of a particular component.  By multiplication by coeff, we can
  64.    adjust the amount of total contribution, giving the formula:
  65.       coeff * (1 - (r/rad)^2)^2
  66.    This varies in strength from coeff at r = 0, to 0 at r = rad. */
  67. void
  68. MakeBlob(blob, threshold, bloblist, npoints, sflag)
  69. BLOB *blob;
  70. DBL threshold;
  71. blobstackptr bloblist;
  72. int npoints;
  73. int sflag;
  74.   {
  75.   unsigned i;
  76.   DBL rad, coeff;
  77.   blobstackptr temp;
  78.   VECTOR mins, maxs;
  79.  
  80.   if (npoints < 1) 
  81.     Error("Need at least one component in a blob.");
  82.  
  83.   blob->threshold = threshold;
  84.   blob->list = (Blob_Element **)malloc(npoints*sizeof(Blob_Element *));
  85.   if (blob->list == NULL) 
  86.     MAError("blob data");
  87.   for (i=0;i < (unsigned)npoints;i++) 
  88.     {
  89.     blob->list[i] = (Blob_Element *)malloc(sizeof(Blob_Element));
  90.     if (blob->list[i] == NULL)
  91.       MAError("blob data");
  92.     }
  93.  
  94.   blob->count = npoints;
  95.   blob->Sturm_Flag = sflag;
  96.  
  97.   /* Initialize the blob data */
  98.   for(i=0;i < (unsigned)npoints;i++) 
  99.     {
  100.     temp = bloblist;
  101.     if (fabs(temp->elem.coeffs[2]) < EPSILON ||
  102.       temp->elem.radius2 < EPSILON) 
  103.       {
  104.       perror("Degenerate blob element\n");
  105.       }
  106.     /* Store blob specific information */
  107.     rad = temp->elem.radius2;
  108.     rad *= rad;
  109.     coeff = temp->elem.coeffs[2];
  110.     blob->list[i]->radius2   = rad;
  111.     blob->list[i]->coeffs[2] = coeff;
  112.     blob->list[i]->coeffs[1] = -(2.0 * coeff) / rad;
  113.     blob->list[i]->coeffs[0] = coeff / (rad * rad);
  114.     blob->list[i]->pos.x = temp->elem.pos.x;
  115.     blob->list[i]->pos.y = temp->elem.pos.y;
  116.     blob->list[i]->pos.z = temp->elem.pos.z;
  117.  
  118.     rad = temp->elem.radius2;
  119.     if (i == 0) 
  120.       {
  121.       /* First component, just set the bounds */
  122.       Make_Vector(&mins,
  123.         temp->elem.pos.x-rad,
  124.         temp->elem.pos.y-rad,
  125.         temp->elem.pos.z-rad)
  126.         Make_Vector(&maxs,
  127.           temp->elem.pos.x+rad,
  128.           temp->elem.pos.y+rad,
  129.           temp->elem.pos.z+rad)
  130.           }
  131.     else 
  132.       {
  133.       /* Check min/max on the bounds */
  134.         mins.x = min(mins.x, temp->elem.pos.x - rad);
  135.       mins.y = min(mins.y, temp->elem.pos.y - rad);
  136.       mins.z = min(mins.z, temp->elem.pos.z - rad);
  137.       maxs.x = max(maxs.x, temp->elem.pos.x + rad);
  138.       maxs.y = max(maxs.y, temp->elem.pos.y + rad);
  139.       maxs.z = max(maxs.z, temp->elem.pos.z + rad);
  140.       }
  141.  
  142.     bloblist = bloblist->next;
  143.     free(temp);
  144.     }
  145.  
  146.   blob->Bounds.Lower_Left = mins;
  147.   VSub(blob->Bounds.Lengths, maxs, mins);
  148.  
  149.   /*  Allocate memory for intersection intervals */
  150.   npoints *= 2;
  151.   blob->intervals = (Blob_Interval *)malloc(npoints*sizeof(Blob_Interval));
  152.   if (blob->intervals == NULL) 
  153.     MAError("blob data");
  154.   }
  155.  
  156. /* Make a sorted list of points along the ray that the various blob
  157.    components start and stop adding their influence.  It would take
  158.    a very complex blob (with many components along the current ray)
  159.    to warrant the overhead of using a faster sort technique. */
  160. static int
  161. determine_influences(P, D, blob, mindist)
  162. VECTOR *P, *D;
  163. BLOB *blob;
  164. DBL mindist;
  165.   {
  166.   int i, j, k, cnt;
  167.   DBL b, t, t0, t1, disc;
  168.   VECTOR V;
  169.   Blob_Interval *intervals = blob->intervals;
  170.  
  171.   cnt = 0;
  172.   for (i=0;i<blob->count;i++) 
  173.     {
  174.     /* Use standard sphere intersection routine
  175.          to determine where the ray hits the volume
  176.          of influence of each component of the blob. */
  177.     VSub(V, blob->list[i]->pos, *P);
  178.     VDot(b, V, *D);
  179.     VDot(t, V, V);
  180.     disc = b * b - t + blob->list[i]->radius2;
  181.     if (disc < EPSILON)
  182.       continue;
  183.     disc = sqrt(disc);
  184.     t1 = b + disc;
  185.     if (t1 < mindist) t1 = 0.0;
  186.     t0 = b - disc;
  187.     if (t0 < mindist) t0 = 0.0;
  188.     if (t1 == t0) continue;
  189.     else if (t1 < t0) 
  190.       {
  191.       disc = t0;
  192.       t0 = t1;
  193.       t1 = disc;
  194.       }
  195.  
  196.     /* Store the points of intersection of this
  197.          blob with the ray.  Keep track of: whether
  198.          this is the start or end point of the hit,
  199.          which component was pierced by the ray,
  200.          and the point along the ray that the
  201.          hit occured at. */
  202.     for (k=0;k<cnt && t0 > intervals[k].bound;k++);
  203.     if (k<cnt) 
  204.       {
  205.       /* This hit point is smaller than one that
  206.             already exists - bump the rest and insert
  207.             it here */
  208.       for (j=cnt;j>k;j--)
  209.         memcpy(&intervals[j], &intervals[j-1],
  210.           sizeof(Blob_Interval));
  211.       intervals[k].type  = 0;
  212.       intervals[k].index = i;
  213.       intervals[k].bound = t0;
  214.       cnt++;
  215.       for (k=k+1;k<cnt && t1 > intervals[k].bound;k++);
  216.       if (k<cnt) 
  217.         {
  218.         for (j=cnt;j>k;j--)
  219.           memcpy(&intervals[j], &intervals[j-1],
  220.             sizeof(Blob_Interval));
  221.         intervals[k].type  = 1;
  222.         intervals[k].index = i;
  223.         intervals[k].bound = t1;
  224.         }
  225.       else 
  226.         {
  227.         intervals[cnt].type  = 1;
  228.         intervals[cnt].index = i;
  229.         intervals[cnt].bound = t1;
  230.         }
  231.       cnt++;
  232.       }
  233.     else 
  234.       {
  235.       /* Just plop the start and end points at
  236.             the end of the list */
  237.       intervals[cnt].type  = 0;
  238.       intervals[cnt].index = i;
  239.       intervals[cnt].bound = t0;
  240.       cnt++;
  241.       intervals[cnt].type  = 1;
  242.       intervals[cnt].index = i;
  243.       intervals[cnt].bound = t1;
  244.       cnt++;
  245.       }
  246.     }
  247.   return cnt;
  248.   }
  249.  
  250.   /* Calculate the field value of a blob - the position vector
  251.    "Pos" must already have been transformed into blob space. */
  252.   static DBL
  253.   calculate_field_value(obj, Pos)
  254.     OBJECT *obj;
  255. VECTOR *Pos;
  256.   {
  257.   int i;
  258.   DBL len, density;
  259.   VECTOR V;
  260.   Blob_Element *ptr;
  261.   BLOB *blob = (BLOB *)obj;
  262.  
  263.   density = 0.0;
  264.   for (i=0;i<blob->count;i++) 
  265.     {
  266.     ptr = blob->list[i];
  267.     VSub(V, ptr->pos, *Pos);
  268.     VDot(len, V, V);
  269.     if (len < ptr->radius2) 
  270.       {
  271.       /* Inside the radius of influence of this
  272.             component, add it's contribution */
  273.       density += len * (len * ptr->coeffs[0] +
  274.         ptr->coeffs[1]) +
  275.       ptr->coeffs[2];
  276.       }
  277.     }
  278.   return density;
  279.   }
  280.  
  281. /* Generate intervals of influence of each component.  After these
  282.    are made, determine their aggregate effect on the ray.  As the
  283.    individual intervals are checked, a quartic is generated
  284.    that represents the density at a particular point on the ray.
  285.  
  286.    After making the substitutions in MakeBlob, there is a formula
  287.    for each component that has the form:
  288.    
  289.       c0 * r^4 + c1 * r^2 + c2.
  290.    
  291.    In order to determine the influence on the ray of all of the
  292.    individual components, we start by determining the distance
  293.    from any point on the ray to the specified point.  This can
  294.    be found using the pythagorean theorem, using C as the center
  295.    of this component, P as the start of the ray, and D as the
  296.    direction of travel of the ray:
  297.  
  298.       r^2 = (t * D + P - C) . (t * D + P - C)
  299.  
  300.    we insert this equation for each appearance of r^2 in the
  301.    components' formula, giving:
  302.  
  303.       r^2 = D.D t^2 + 2 t D . (P - C) + (P - C) . (P - C)
  304.  
  305.    Since the direction vector has been normalized, D.D = 1.
  306.    Using the substitutions:
  307.  
  308.       t0 = (P - C) . (P - C),
  309.       t1 = D . (P - C)
  310.  
  311.    We can write the formula as:
  312.  
  313.       r^2 = t0 + 2 t t1 + t^2
  314.  
  315.    Taking r^2 and substituting into the formula for this component
  316.    of the blob we get the formula:
  317.  
  318.       density = c0 * (r^2)^2 + c1 * r^2 + c2,
  319.  
  320.    or:
  321.  
  322.       density = c0 * (t0 + 2 t t1 + t^2)^2 +
  323.                 c1 * (t0 + 2 t t1 + t^2) +
  324.                 c2
  325.  
  326.    Expanding terms and collecting with respect to "t" gives:
  327.       t^4 * c0 +
  328.       t^3 * 4 c0 t1 +
  329.       t^2 * (c1 + 2 * c0 t0 + 4 c0 t1^2)
  330.       t   * 2 (c1 t1 + 2 c0 t0 t1) +
  331.             c2 + c1*t0 + c0*t0^2
  332.  
  333.    This formula can now be solved for "t" by any of the quartic
  334.    root solvers that are available.
  335. */
  336. int All_Blob_Intersections(Object, Ray, Depth_Stack)
  337. OBJECT *Object;
  338. RAY *Ray;
  339. ISTACK *Depth_Stack;
  340.   {
  341.   BLOB *blob = (BLOB *)Object;
  342.   DBL dist, len, *tcoeffs, coeffs[5], roots[4];
  343.   int i, j, cnt;
  344.   VECTOR P, D, V;
  345.   int root_count, in_flag;
  346.   Blob_Element *element;
  347.   DBL t0, t1, c0, c1, c2;
  348.   VECTOR IPoint, dv;
  349.   Blob_Interval *intervals = blob->intervals;
  350.   int Intersection_Found = FALSE;
  351.  
  352.   Ray_Blob_Tests++;
  353.  
  354.   /* Transform the ray into the blob space */
  355.   if (blob->Trans != NULL) 
  356.     {
  357.     MInvTransPoint(&P, &Ray->Initial, blob->Trans);
  358.     MInvTransDirection(&D, &Ray->Direction, blob->Trans);
  359.     }
  360.   else 
  361.     {
  362.     P = Ray->Initial;
  363.     D = Ray->Direction;
  364.     }
  365.  
  366.     len = sqrt(D.x * D.x + D.y * D.y + D.z * D.z);
  367.   if (len == 0.0)
  368.     return 0;
  369.   else 
  370.     {
  371.     D.x /= len;
  372.     D.y /= len;
  373.     D.z /= len;
  374.     }
  375.  
  376.     /* Figure out the intervals along the ray where each
  377.    component of the blob has an effect. */
  378.     if ((cnt = determine_influences(&P, &D, blob, 0.01)) == 0)
  379.       /* Ray doesn't hit the sphere of influence of any of
  380.       its component elements */
  381.       return 0;
  382.  
  383.   /* Clear out the coefficients */
  384.   for (i=0;i<4;i++) coeffs[i] = 0.0;
  385.   coeffs[4] = -blob->threshold;
  386.  
  387.   /* Step through the list of influence points, adding the
  388.       influence of each blob component as it appears */
  389.   for (i=0,in_flag=0;i<cnt;i++) 
  390.     {
  391.     if (intervals[i].type == 0) 
  392.       {
  393.       /* Something is just starting to influence the ray,
  394.             so calculate its coefficients and add them
  395.          into the pot. */
  396.       in_flag++;
  397.       element = blob->list[intervals[i].index];
  398.  
  399.       VSub(V, P, element->pos);
  400.       c0 = element->coeffs[0];
  401.       c1 = element->coeffs[1];
  402.       c2 = element->coeffs[2];
  403.       VDot(t0, V, V);
  404.       VDot(t1, V, D);
  405.       tcoeffs = &(element->tcoeffs[0]);
  406.  
  407.       tcoeffs[0] = c0;
  408.       tcoeffs[1] = 4.0 * c0 * t1;
  409.       tcoeffs[2] = 2.0 * c0 * (2.0 * t1 * t1 + t0) + c1;
  410.       tcoeffs[3] = 2.0 * t1 * (2.0 * c0 * t0 + c1);
  411.       tcoeffs[4] = c0 * t0 * t0 + c1 * t0 + c2;
  412.  
  413.       for (j=0;j<5;j++) coeffs[j] += tcoeffs[j];
  414.       }
  415.     else 
  416.       {
  417.       /* We are losing the influence of a component, so
  418.             subtract off its coefficients */
  419.       tcoeffs = &(blob->list[intervals[i].index]->tcoeffs[0]);
  420.       for (j=0;j<5;j++) coeffs[j] -= tcoeffs[j];
  421.       if (--in_flag == 0)
  422.         /* None of the components are currently affecting
  423.                the ray - skip ahead. */
  424.         continue;
  425.       }
  426.  
  427.     /* Figure out which root solver to use */
  428.     if (blob->Sturm_Flag == 0)
  429.       /* Use Ferrari's method */
  430.       root_count = solve_quartic(coeffs, &roots[0]);
  431.     else
  432.       /* Sturm sequences */
  433.       if (fabs(coeffs[0]) < COEFF_LIMIT)
  434.         if (fabs(coeffs[1]) < COEFF_LIMIT)
  435.           root_count = solve_quadratic(&coeffs[2], &roots[0]);
  436.         else
  437.           root_count = polysolve(3, &coeffs[1], &roots[0]);
  438.       else
  439.         root_count = polysolve(4, coeffs, &roots[0]);
  440.  
  441.     /* See if any of the roots are valid */
  442.     for(j=0;j<root_count;j++) 
  443.       {
  444.       dist = roots[j];
  445.       /* First see if the root is in the interval of influence of
  446.             the currently active components of the blob */
  447.       if ((dist >= intervals[i].bound) &&
  448.         (dist <= intervals[i+1].bound) &&
  449.         (dist > Blob_Tolerance)) 
  450.         {
  451.         VScale(IPoint, D, dist);
  452.         VAdd(IPoint, IPoint, P);
  453.         /* Transform the point into world space */
  454.         if (blob->Trans != NULL)
  455.           MTransPoint(&IPoint, &IPoint, blob->Trans);
  456.         VSub(dv, IPoint, Ray->Initial);
  457.         VLength(len, dv);
  458.         if (Point_In_Clip(&IPoint, Object->Clip)) 
  459.           {
  460.           push_entry(len,IPoint,Object,Depth_Stack);
  461.           Intersection_Found = TRUE;
  462.           }
  463.         }
  464.       }
  465.     }
  466.   if(Intersection_Found)
  467.     Ray_Blob_Tests_Succeeded++;
  468.   return Intersection_Found;
  469.   }
  470.  
  471. /* Calculate the density at this point, then compare to
  472.    the threshold to see if we are in or out of the blob */
  473. int Inside_Blob (Test_Point, Object)
  474. VECTOR *Test_Point;
  475. OBJECT *Object;
  476.   {
  477.   VECTOR New_Point;
  478.   BLOB *blob = (BLOB *) Object;
  479.  
  480.   /* Transform the point into blob space */
  481.   if (blob->Trans != NULL)
  482.     MInvTransPoint(&New_Point, Test_Point, blob->Trans);
  483.   else
  484.     New_Point = *Test_Point;
  485.  
  486.   if (calculate_field_value(Object, &New_Point) >
  487.     blob->threshold-INSIDE_TOLERANCE)
  488.     return ((int) 1-blob->Inverted);
  489.   else
  490.     return ((int) blob->Inverted);
  491.   }
  492.  
  493. void Blob_Normal (Result, Object, IPoint)
  494. OBJECT *Object;
  495. VECTOR *Result, *IPoint;
  496.   {
  497.   VECTOR New_Point, V;
  498.   int i;
  499.   DBL dist, val;
  500.   BLOB *blob = (BLOB *) Object;
  501.   Blob_Element *temp;
  502.  
  503.   /* Transform the point into the blobs space */
  504.   if (blob->Trans != NULL)
  505.     MInvTransPoint(&New_Point, IPoint, blob->Trans);
  506.   else
  507.     New_Point = *IPoint;
  508.  
  509.   Make_Vector(Result, 0, 0, 0);
  510.  
  511.   /* For each component that contributes to this point, add
  512.       its bit to the normal */
  513.   for(i=0;i<blob->count;i++) 
  514.     {
  515.     temp = blob->list[i];
  516.     V.x = New_Point.x - temp->pos.x;
  517.     V.y = New_Point.y - temp->pos.y;
  518.     V.z = New_Point.z - temp->pos.z;
  519.     dist = (V.x * V.x + V.y * V.y + V.z * V.z);
  520.  
  521.     if (dist <= temp->radius2) 
  522.       {
  523.       val = -2.0 * (2.0 * temp->coeffs[0] * dist +
  524.         temp->coeffs[1]);
  525.       Result->x += val * V.x;
  526.       Result->y += val * V.y;
  527.       Result->z += val * V.z;
  528.       }
  529.     }
  530.   val = (Result->x * Result->x + Result->y * Result->y +
  531.     Result->z * Result->z);
  532.   if (val < EPSILON) 
  533.     {
  534.     Result->x = 1.0;
  535.     Result->y = 0.0;
  536.     Result->z = 0.0;
  537.     }
  538.   else 
  539.     {
  540.     val = 1.0 / sqrt(val);
  541.     Result->x *= val;
  542.     Result->y *= val;
  543.     Result->z *= val;
  544.     }
  545.  
  546.     /* Transform back to world space */
  547.     if (blob->Trans != NULL)
  548.       MTransNormal(Result, Result, blob->Trans);
  549.   VNormalize(*Result, *Result);
  550.   }
  551.  
  552. void Translate_Blob (Object, Vector)
  553. OBJECT *Object;
  554. VECTOR *Vector;
  555.   {
  556.   TRANSFORM Trans;
  557.   Compute_Translation_Transform(&Trans, Vector);
  558.   Transform_Blob(Object, &Trans);
  559.   }
  560.  
  561. void Rotate_Blob (Object, Vector)
  562. OBJECT *Object;
  563. VECTOR *Vector;
  564.   {
  565.   TRANSFORM Trans;
  566.   Compute_Rotation_Transform(&Trans, Vector);
  567.   Transform_Blob(Object, &Trans);
  568.   }
  569.  
  570. void Scale_Blob (Object, Vector)
  571. OBJECT *Object;
  572. VECTOR *Vector;
  573.   {
  574.   TRANSFORM Trans;
  575.   Compute_Scaling_Transform(&Trans, Vector);
  576.   Transform_Blob(Object, &Trans);
  577.   }
  578.  
  579. void Transform_Blob(Object,Trans)
  580. OBJECT *Object;
  581. TRANSFORM *Trans;
  582.   {
  583.   if (((BLOB *)Object)->Trans == NULL)
  584.     ((BLOB *)Object)->Trans = Create_Transform();
  585.   recompute_bbox(&Object->Bounds, Trans);
  586.   Compose_Transforms(((BLOB *)Object)->Trans, Trans);
  587.   }
  588.  
  589. void Invert_Blob(Object)
  590. OBJECT *Object;
  591.   {
  592.   ((BLOB *) Object)->Inverted = 1 - ((BLOB *)Object)->Inverted;
  593.   }
  594.  
  595. void *Copy_Blob (Object)
  596. OBJECT *Object;
  597.   {
  598.   BLOB *New;
  599.   unsigned i, cnt;
  600.  
  601.   New = Create_Blob();
  602.   *New = * ((BLOB *)Object);
  603.  
  604.   New->Trans = Copy_Transform(New->Trans);
  605.  
  606.   cnt = ((BLOB *)Object)->count;
  607.   New->list = (Blob_Element **)malloc(cnt * sizeof(Blob_Element *));
  608.   if (New->list == NULL) 
  609.     MAError("blob data");
  610.   for (i=0;i<cnt;i++) 
  611.     {
  612.     New->list[i] = (Blob_Element *)malloc(sizeof(Blob_Element));
  613.     if (New->list[i] == NULL)
  614.       MAError("blob data");
  615.     memcpy(New->list[i], ((BLOB *)Object)->list[i], sizeof(Blob_Element));
  616.     }
  617.  
  618.   New->intervals = (Blob_Interval *)malloc(2*New->count*sizeof(Blob_Interval));
  619.   if (New->intervals == NULL)
  620.     MAError("blob data");
  621.  
  622.   return (New);
  623.   }
  624.  
  625. /* Allocate a blob. */
  626. BLOB *Create_Blob()
  627.   {
  628.   BLOB *New;
  629.  
  630.   if ((New = (BLOB *) malloc (sizeof (BLOB))) == NULL)
  631.     MAError ("blob");
  632.  
  633.   INIT_OBJECT_FIELDS(New, BLOB_OBJECT, &Blob_Methods)
  634.  
  635.     New->Trans = NULL;
  636.   New->Inverted = FALSE;
  637.   New->count = 0;
  638.   New->threshold = 0.0;
  639.   New->list = NULL;
  640.   New->intervals = NULL;
  641.   New->Sturm_Flag = FALSE;
  642.  
  643.   return (New);
  644.   }
  645.  
  646. void Destroy_Blob (Object)
  647. OBJECT *Object;
  648.   {
  649.   unsigned i;
  650.  
  651.   Destroy_Transform(((BLOB *)Object)->Trans);
  652.   for (i=0;i < (unsigned)((BLOB *)Object)->count;i++)
  653.     free (((BLOB *)Object)->list[i]);
  654.   free (((BLOB *)Object)->list);
  655.   free (((BLOB *)Object)->intervals);
  656.   free (Object);
  657.   }
  658.